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SUMMARY 

A parameter and structure estimation technique for multivariable systems is used 
to obtain a state space representation of open loop dynamics of the space shuttle 
main engine in a-canonical form. The parametrizat ion being used is both minimal and 
unique. The simplified linear models may be used for fault detection studies and 
control system design and development. 


INTRODUCTION 

Accurate representation of the dynamic behavior of the space shuttle main engine 
( SSME ) is required for a variety of diagnostic, control and evaluation purposes. A 
complete nonlinear dynamic simulation has been developed (Rocketdyne Division of 
Rockwell International Corporation, 1981). However, its size and complexity make it 
very difficult to use for design purposes. In a previous study, Duyar, Guo and 
Merrill, 1990, used a single-input single-output (SISO) identification technique to 
obtain a transfer function representation of the dynamic behavior of the SSME at 
100 percent power level. Although SISO identification approach is quite simple, it 
results in a nonminimal representation of the SSME. This paper extends the previous 
work by developing linearized dynamic models of the SSME at five different power 
levels. A multivariable identification and a minimal realization technique (Eldem 
and Yildizbayrak, 1988; Eldem and Duyar, 1989) is used to present these models in 
a-canonical form. The paper also develops an identification technique for piecewise 
linear systems with static nonlinear gains, which is used in modeling the dynamics of 
the SSME. 

The identified linearized models are valid in a limited response region about 
the operating points corresponding to the power levels. The models are useful for 
real time estimation and fault detection as well as open loop engine dynamic studies 
and closed loop control analysis using a user generated control law. 

Initially a brief description of the SSME is given. This is followed by a 
description of the identification scheme and the model used. Finally, results 
obtained from the identified models are compared with the results obtained form the 
nonlinear simulation for the same input. 


NOMENCLATURE 


A,B,C, 

k 

MR 

N 

P 

S 

T 

t 

t 


state space matrices 
discrete time integer 

oxygen to fuel mixture ratio of the main combustion chamber 

number of samples 

pressure 

speed 

temperature 

time 

clock time 


Greek Symbols: 

P valve actuator output rotary motion 

<5u deviation of input from nominal 

<5x deviation of state from nominal 

6y deviation of output from nominal 

X gain 

Subscripts : 

C chamber 

CCV coolant control valve 

FPOV fuel preburner oxidizer valve 

HPFT high pressure fuel turbine 

HPOT high pressure oxidizer turbine 

MFV main fuel valve 

MOV main oxidizer valve 

OPOV oxidizer preburner oxidizer valve 


THE SPACE SHUTTLE MAIN ENGINE 

Different aspects of the SSME as well as its principles of operation are 
described in the literature (Duyar, Guo and Merrill, 1990; Landauer, 1988; Klatt 
et al., 1982). For the sake of completeness a brief description of the main engine 
is also given below. 

The space shuttle orbiter main propulsion system is composed of three main 
engines. The engines use liquid oxygen and liquid hydrogen propellants carried in an 
external tank attached to the orbiter. To understand the overall flow of fuel and 
oxidizer to produce the thrust, a schematic diagram of the propellant flows and the 
control valves is shown in figure 1. 

The two high pressure turbines are driven by a fuel turbine preburner and an 
oxidizer turbine preburner, each of which produces hot gas. The low pressure 
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turbines are driven by the high pressure pump flows. The fuel from the high pressure 
fuel pump ( HPFP ) goes through the main fuel valve (MFV) . After the MFV, the flow 
divides into fixed nozzle cooling flow, main chamber cooling flow, and coolant 
control valve (CCV) flow. Heat is absorbed from the combustion chamber and nozzle. 
The fixed nozzle cooling and the coolant control valve flows then recombine and 
travel to the preburners where combustion and pressure are controlled by the fuel 
preburner oxidizer valve (FPOV) and the oxidizer preburner oxidizer valve (OPOV) . 

Main chamber coolant flow is used to drive the low pressure fuel turbopump. 


The existing control system of the SSME uses five valves; FPOV, OPOV, MFV, MOV, 
and CCV. They control the mixture ratio and the main chamber pressure which is 
correlated with the thrust. Open and closed loop control of these valves are used to 
accomplish the SSME mission. A typical SSME mission is shown in figure 2. For the 
purpose of model identification of the open loop system, the rotary motion of the 
five valves (/? opov , /? Fpov , /? ccv , /? M0V , is used as input. The outputs are the 

chamber pressure, P c , and the mixture ratio, MR. 


In a previous study, Duyar, Guo and Merrill, 1990, showed that (3 , , and 

P valve rotary motions are essentially decoupled from the outputs during the main 


stage of operation. In this study, only the valve rotary motions, p 
are used as the inputs of the open loop system. 


FPOV 


and B 

r ( 


OPOV' 


SYSTEM IDENTIFICATION 

According to the classical definition given by Zadeh, 1962, identification is 
the determination, on the basis of input and output, of a system within a specified 
class of systems, to which the system under test is equivalent. This definition is 
quite general and allows many degrees of freedom in the practical formulation of the 
identification problem. For example, one may select a linear or a nonlinear model. 

It is possible to have different parametrizations for a given model. Additional 
choices need to be made to select the input signal. In this section a comparative 
review of the literature is given to make these choices. 

There are two major classes of mathematical models which are used for 
identification: linear and nonlinear models. Due to the complexity inherent in 

nonlinear systems and due to the fact that there is no well established general 
theory for such systems, it is assumed in this study that the physical process being 
considered behaves linearly at least within a small neighborhood of a nominal 
operating point. 

Parametrization is an important issue for multi-input multi-output (MIMO) 
systems since they admit more than one parametrization depending on their 
observability indices. Thus, in order to obtain a minimal parametrization for MIMO 
systems, the structure of the system i.e., the observability indices related to each 
output must be determined. The detailed analysis of structure identification can be 
found in the literature (Guidorzi, 1975, 1981; Wertz et al., 1982; Overbeek and 
Ljung, 1982; Beghelli and Guidorzi, 1983; and Correa and Glover, 1984). Recently, 
Eldem and Yildizbayrak, 1988, employed the notion of output injection to obtain a new 
canonical form for a special class of observable systems. Through this canonical 
form, they developed a structure and parameter identification algorithm for a 
restricted class of systems. Eldem and Duyar, 1989, extended this technique to the 
class of all minimal systems and constrained the output injection by an extra 
condition to guarantee the uniqueness of the parametrization. In this study, the 
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techniques developed by Eldem and Duyar, 1989, are adopted for the identification of 
the space shuttle main engine dynamics. 

Selection of an appropriate input signal is an important step in identification 
problems. A basic criterion for this selection is that the input/output data should 
be informative enough to discriminate between different models among the class of 
models being considered (Ljung, 1987). Without this discrimination there is no 
guarantee that the obtained parameters are the true parameters of the system. This 
criterion can be expressed mathematically in terms of the covariance matrices of the 
input signals and the order of the system being identified (Ljung, 1987; Bitmead, 
1984; and Anderson and Johnson, 1982). Uncorrelated pseudo random binary sequences 
which are used in this study are examples of such input signals which can be used for 
designing an informative experiment. 

Based on the above review, the steps in identification consist of: 

(1) selection of a driving signal with persistent excitations, (2) selection of a 
model, (3) parameter and structure estimation, and (4) model verification. These 
steps are followed for the identification of the SSME and are outlined below. 


Selection of a Driving Signal 

Pseudo random binary sequences (PRBS) are selected as the input perturbation 
signals to excite the system because of their convenience and suitability for similar 
applications (Cottington and Pease, 1978; Sudhakar et al., 1988). Two uncorrelated 
sequences, which satisfy the requirement of persistent excitation, are used for OPOV 
and FPOV inputs. The sequences have a clock time of 0.04 sec and a length of 127. 
This corresponds to a maximum frequency on 78.5 rad/sec (12.5 Hz), a minimum 
frequency of 1.24 rad/sec (0.2 Hz) and a signal duration of 5.08 sec. 


Selection of a Mathematical Model 

In a previous study by Duyar, Guo and Merrill, 1990, the existence of 
significant backlash and stiction nonlinearities associated with the valve dynamics 
were observed. These two nonlinearities can be isolated from the rest of the system 
and dealt with separately. Here, the open loop system includes actuator dynamics, 
the nonlinear element which contains valve stiction and valve linkage backlash, and 
the engine system which includes all other nonlinearities and the engine dynamics. 

To simplify the analysis, the backlash and stiction are temporarily removed from the 
analysis. Now the open loop system, without the backlash and the stiction, is 
linearized and identified. The backlash and the stiction nonlinearities can later be 
added to the identified system. In practice values for backlash and stiction are 
determined from manufacturer's specifications or bench testing of the individual 
components . 

The nonlinear dynamics, except the backlash and the stiction, of the SSME can be 
described by the nonlinear equations: 


x(t) = f [x(t) , u(t) ] 


where x, u, and y are the state, the control and the output vectors, respectively. 
Linearizing these equations about a nominal operating condition and discretizing 
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y (t ) = g[x(t) ] 


( 2 ) 


yields : 

<5x(k + 1) = A <5x ( k) + B <5u ( k ) ( 3 ) 

<5y ( k) = C <5x(k) ( 4 ) 

where <5x, <5u, and 5 y are the deviations of the state, the input, and the output 
vectors about the nominal operating condition. It is assumed that the system 
described by equations (3) and (4) is stable and observable and the C matrix has 
full row rank. Furthermore, this system is realized in a-canonical form (Eldem and 
Yildizbayrak, 1988; and Eldem and Duyar, 1989), i.e., the following relations hold: 

C = [o s H -1 ] < 5 > 


A = A + KHC , 

o 


with A^ = 0 

o 


( 6 ) 


(HC) 



= 0 


(7) 


( HC )r. A o K c. = °' for /*i > /*j and k < Mi - (8) 

Here the subscripts and c^ denote the i'th row and j'th column 

respectively. Superscripts indicate exponentiation. The structure matrix A is 
lower left triangular and consists of zeros and ones only and is determined by the 
observability indices, ji i where i associates with the i'th output and 

fi = max The matrix K is a deadbeat observer gain. 

In the following section a parameter and structure identification scheme for the 
above system is given. It is assumed that H matrix is equal the identity matrix, 

I. The generalization of this identification scheme to the case with H t I is 
treated by Eldem and Duyar, 1989. 


Identification of Linear Models 

Throughout this section a linear, discrete time, multivariable system described 
by equations (3) and (4) is considered. The triplet (C,A,B) satisfies equations (5) 
to (8) with H = I. The expression for the state <5x(k) at time k of this system 
is given as: 
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<5x(k) = h o 6x(0) + £ \ -1 [K B] 


i=l 


<5y(k - i) 
<5u(k - i) 


(9) 


Using the nilpotency of A the above equation yields: 


5x(k) = £ A^[K B] 


i=l 


<5y(k - i) 
<5u(k - i) 


, for k > H 


( 10 ) 


This implies that 


i-l 


5y(k) = Y, CA o t K B 1 


i=l 


<5y(k - i) 
6u ( k - i) 


, for k > fi 


( 11 ) 


Now let 


= 


CK CA K 


M- 1 

CA K 


( 12 ) 


6 u = 


CB CA B . 

o 


fj.-l 

. CA B 

o 


(13) 


<Py - [<5y(k - 1) T ... 5y(k - /i ) T ] T 


(14) 


K = [<5 u ( k - 1) T ... <5u(k - /i) T ] T 


(15) 


and 


0 = \e y e a 


(16) 


0(k) = ^y(k) T (Mk) T 


(17) 


Using the notation above , equation (11) can be written as: 

<5y(k) = 0<j){ k) , for k > fl 

If the length of the input/output data is N, then: 

Y n = [<5y(/i) 6y(/i +1) ... <5y(N) ] 


(18) 


(19) 
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</> Y (H) ... <f) Y ( N) 

<P a (H) ... <f>JU) 


( 20 ) 


Defining Q and R as 


Q = ...0 y (N)] 


( 21 ) 


R = [^„(/0 ^ U ( N )] 


(22) 


we obtain 


Y N = do + e u R 


(23) 


and postmultiplying the above equation with [q t R t ] yields: 

= y n[q T rT ] 


(24) 


7 n = K e 




T T 

QQ QR 


T T 

RQ RR 


(25) 


Defining $ as: 


$ = 


T T 

QQ QR 


T T 

RQ RR 


(26) 


Equation (25) can be rewritten as: 


y n = 0$ 


(27) 


Equation (27) can be solved to determine the parameters 6 and 6 if the 
observability indices /i. are known. To determine them, an expression for the 
observability matrix is derived. Toward this end first a matrix Q q is defined as: 


Q 0 = QQ T - QR T (rR T ) _ 1 RQ T 


(28) 


Under the assumption that inputs are "persistently exciting" RR T is nonsingular. 
Also let V and U n _ +i be the observability and controllability matrices written 
for ji and^ N-/J + 1 s£eps respectively, i.e.. 
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v M = [ca m ~ 1 ca m ' 2 . . . c] t 


= 


B AB . . . A 




(29] 


(30) 


Then it has been shown by Eldem and Duyar (1989) that the following expression holds: 

Qo=y< < 31 > 


Here M is a nonsingular matrix defined as: 


M 


U N-/i+ l R o 



-l 


T T 

R U 
o 


(32) 


and R q is given below: 


R o = 


<5u(0) <5u(l) 
o 5u(0) 


5u(n - n) 
5u(n - n - l) 

5u(0) 


(33) 


Equation (31) has an important consequence: row by row rank search on Q q yields 

the observability indices of the system. This is due to the fact that post- 
multiplying the observability matrix V M by MV T does not change the observability 

T A* 

indices since M is nonsingular and V has full row rank. Thus equation (27) can 
now be solved for the parameter matrix 0, using the least square technique. 


Identification of Piecewise Linear Models 

The procedure outlined earlier is used to identify the open loop system dynamics 
of the SSME at different power levels. The responses of the identified open loop 
system are compared with the responses obtained from the nonlinear simulation. In 
general , it is observed that the gains for positive and negative perturbation signals 
are significantly different. In order to compensate for this phenomenon, a system 
gain with different values for positive and negative perturbations is added to the 
linearized model. 

With this modification the system equations become 

<5x(k + 1) = A <5x(k) + B <5u*(k) (34) 


6y(k) = C 5x(k) 


(35) 
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where <5u (k) is defined component-wise as 



X i 5u i (k) , for 5u i (k) > 0 

(2 - X i )6u.(k) , for Su^k) < 0 


The input sequence <5u(k) can now be decomposed as 

6u(k) = <5u + (k) + <5u (k) 


where 



<5u.(k) , for 5u i (k) > 0 
0 , otherwise 


(36) 


(37) 


(38) 


5u i (k) 


<5u.(k) , for <5u.(k) < 0 
0 , otherwise 


(39) 


Using the above definitions the system equations can be written as 

<5x(k + 1) = A <5x(k) + B 5u(k) + B(A - I)[<5u + (k) - <5u ( k ) ] 
<5y ( k) = C <5x(k) 


(40) 


where A = diag {X.}. Using a-canonical form (A = A° + KC) and the properties of 
A q an expression for <5y(k) can be written as 


<5y(k) = 0 y <p Y [ k) + G u <f> u { k) 


+ 0 U (A - i) 


0 u < k > " ^u< k > 


(41) 


For constant amplitude perturbations (for instance PRBS) the last term in the above 
equation is constant since 


6wr(k) - (5u 1 (k) = | 5u i (k) | 


Thus, taking the ensemble averages in equation (41) yields 

6 *. = 0y 5 y m + ^U 5u ,, + ^ U (A - I) 7 e 


(42) 
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where 5y m and 6 u^ denote the average values of u(*) and y(*). 
Furthermore, 


respectively. 


9 = V C A 1 ” 1 KH 

y o 

i=l 


( 43 ) 


M 

* u = £ 


i= 1 


i-1 

C A B 

o 


(44) 


e is an m-vector consisting of ones and 7 is the magnitude of the perturbations. 
Solving equation (42) for X.'s we obtain 


X 





<5u 


(45) 


where X = X 2 ... , e = [ll...l] T . 

Note that using equation (42) and the fact that u + (k) - u’(k) = constant, the 
piecewise linear system defined by equations (34) to (36) can be transformed into a 
purely linear system as follows: 


Let x(k) = <5x(k) - Sx^ 
y ( k) = 5y(k) - (5y m 
u(k) = <5u(k) - <5u m 


Then, by equations (40) and (42), it follows that 

x ( k + 1) = Ax ( k) + Bu(k) < 46 ) 

y (k) = Cx ( k) < 47 ) 

Now the identification procedure outlined in this section can be used to 
determine (C,A,B) of the above system. Since C,A,B is assumed to be in a- 

canonical form, i.e., A = A q + KC, and 9 U can be calculated by equations (43) 

and (44). Then, using the ensemble averages 5 y^ and 6 u^ in equation (45) we can 
obtain the system gains {X i >. 


RESULTS 

Using the above procedure and the identification algorithm presented in this 
paper, the parameters of the piecewise linear model described by equations ( 34 ) 
to (36) are identified at 70, 80, 90, 100, and 110 percent power levels. As 
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mentioned earlier, the input vector, <5u, represents the deviations of the valve 
rotary motion from the nominal values defined as: 


«“ ■ [%,ov f 


The output vector <5 y represents the deviations of the chamber inlet pressure and 
the mixture ratio from their nominal values and defined as: 

<5y = [<5 p c 6mr] t < 49 > 


The A, B, C matrices, the nonlinear gains X i and X 2 and the observability indices 
of the identified piecewise linear model are given in table I. 


The validity of the estimated parameters of the system is checked by comparing 
the responses obtained from the identified system with the response of the nonlinear 

simulation. Both a state variable filter and the model of the identified system are 
used for comparison purposes. The state variable filter utilizes the measurements of 
both the output <5 y, and the input, <5u, to estimate the next value of the output, 
<5y f (k), defined as 




5 y t = E 

i=l 


i-1 

CA [KB] 
o L J 


5y(k - i) 
5u*( k - i) 


( 50 ) 


The model of the identified system utilizes only the measurement of the input data, 
<5u, to predict the output, 5 y , and the state, <5x, as 


<5x ( k + 1) = A<5 x ( k) + B<5u*(k) 


( 51 ) 


<5y m (k) = c<5x(k) 


( 52 ) 


Two different test input signals are used for comparison purposes. The first 
test input signal consists of two full length PRBS (different PRBS than the ones used 
for identification purposes) which cover the same frequency range as those used for 
identification. The second test signal consists of step and ramp inputs as shown in 
figure 3, and covers a lower frequency range than the range of validity of the 
identified system. 

In order to compare the responses of the identified system with the responses 
obtained from the nonlinear simulation, a standard error of estimate (SEE) is defined 
as 
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SEE 


[5z.(k) - <5y.(k) 


\ 


E 5y i (k)‘ 

k= 1 


( 53 ) 


Here the subscript i denotes the i th element of the output vector <5 y and 6 z 
refers either to the model output, i.e., 5 z = 5 y or the filter output, i.e., 

<5z = 5 y f . The SEE's obtained for both of the test signals, given in table II, 
indicate very good agreement. The comparison of the responses of both the identified 
model and the filter to the responses obtained from the nonlinear simulation of the 
SSME for low frequency test signals also indicate good agreement as shown in 
figures 4 to 8. 


CONCLUSION 

A multivariable system identification technique is used to represent the dynamic 
behavior of the SSME at five different power levels. Minimal realizations of these 
models in state space representation using a-canonical form are given. An 
identification technique to determine static nonlinear gain changes is also 
presented . 

The comparison of the responses of the nonlinear simulation with the responses 
of identified models indicates very good agreement and can be used for control design 
purposes. Since the identified models are valid for limited response regions, this 
study will be extended so as to obtain models at other operating points. Then these 
point models will be linked to obtain a simplified model of the SSME covering its 
full range of operation. 
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TABLE I. - IDENTIFIED SYSTEM PARAMETERS AT 70, 80, 90, 100, and 110 PERCENT POWER LEVELS 


Power 
level , 
percent 

A 

B 

C 

K 

V 

70 

0 

0 

-0.2553 

0.0416 

0.3923 

0.0958 








0 

0 

0.5044 

-0.1108 

1.0558 

-0.2652 

0 

0 

1 

0 

0.8957 

0.8403 


1 

0 

1.0218 

-0.113 

0.2584 

0.0788 

0 

0 

0 

1 




0 

1 

-0.7902 

0.5812 

0.5902 

-0.2697 







80 

0 

0 

-0.1576 

0.0257 

0.4068 

0.1564 








0 

0 

0.2587 

-0.0545 

0.9631 

-0.674 

0 

0 

1 

0 

1.1211 

0.9279 


1 

0 

0.8759 

-0.0666 

0.2574 

0.0384 

0 

0 

0 

1 




0 

1 

-0.5157 

0.4766 

0.6046 

-0.2307 







90 

0 

0 

-0.1206 

0.0155 

0.4223 

0.1546 








0 

0 

0.1434 

-0.0275 

0.9217 

-0.7302 

0 

0 

1 

0 

0.9386 

0.7699 


1 

0 

0.8117 

-0.0478 

0.2518 

0.021 

0 

0 

0 

1 




0 

1 

-0.4183 

0.4226 

0.6288 

-0.1779 







100 

0 

0 

-0.1279 

0.0212 

0.3534 

0.0475 








0 

0 

0.1994 

-0.0508 

0.6831 

-0.6446 

0 

0 

1 

0 

0.8594 

0.9047 


1 

0 

0.7816 

-0.0565 

0.2007 

-0.0181 

0 

0 

0 

1 




0 

1 

-0.5253 

0.475 

0.5043 

-0.1387 







110 

0 

0 

-0.084 

0.0026 

0.2416 

0.0786 








0 

0 

0.1082 

-0.0123 

0.5112 

-0.575 

0 

0 

1 

0 

0.9704 

1.1231 


1 

0 

0.728 

0.0181 

0.1456 

-0.0243 

0 

0 

0 

1 




0 

1 

-0.4164 

0.333 

0.3593 

-0.078 








TABLE II. - STANDARD ERROR OF ESTIMATES (SEE) OF 


PREDICTIONS OF THE FILTER AND THE MODEL 



Operating 

condition 



PRBS signal 


Power, percent 


70 

80 

90 

100 

110 

Filter 

P 

c 

MR 

Model 

P 

c 

MR 

0.0648 

0.0876 

0.0793 

0.0667 

0.0488 

0.1028 

0.0951 

0.0943 

0.0443 
0 . 088 

0.0819 

0.0816 

0.0441 

0.0526 

0.0781 

0.0538 

0.0985 

0.0822 

0.124 

0.0903 


Low frequency test signal 


Power, percent 


70 

80 

90 

100 

110 

Filter 

P c 

MR 

Model 

P c 

MR 

0.0123 

0.0483 

0.0481 

0.0911 

0.0422 

0.0721 

0.14 

0.0858 

0.0508 

0.1024 

0.1536 

0.1337 

0.0356 

0.074 

0.0926 

0.1002 

0.0204 

0.1049 

0.0469 

0.1639 
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Power level in percent, RPL 
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Figure 1 —Propellant flow and control valves (Rockwell, 1988). 




0 2 4 6 8 10 12 

Time, sec 

Figure 3. — Low frequency test input signals, 8 OPO v 
and B p PO v 
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Normalized output 


.06 




(b) MR, mixture ratio 

Figure 4. — Comparison of the responses of the identified model and the state variable 
filter with the nonlinear simulation a 70% power lever. 
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Normalized output 




(b) MR mixture ratio. 


Figure 5. — Comparison of the responses of the identified model and the state variable filter 
with the nonlinear simulation a 80% power lever. 
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Normalized output 



(a) P c , chamber inlet pressure. 



(b) MR mixture ratio. 

Figure 6. — Comparison of the responses of the identified model and the state variable filter 
with the nonlinear simulation a 90% power lever. 
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Normalized output 


.04 


Simulation 

Filter 

Model 




Time, sec 

(b) MR, mixture ratio 


Figure 7. — Comparison of the responses of the identified model and the state variable 
filter with the nonlinear simulation a 1 00% power lever. 
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Normalized output 


.03 



-.03 


(a) P c , chamber inlet pressure. 


J i 



(b) MR, mixture ratio 

Figure 8. — Comparison of the responses of the identified model and the state variable 
filter with the nonlinear simulation a 1 1 0% power lever. 
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